Aggregation of self-propelled colloidal rods near confining walls 
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Non-equilibrium collective behavior of self-propelled colloidal rods in a confining channel is studied 
using Brownian dynamics simulations and dynamical density functional theory. We observe an 
aggregation process in which rods self-organize into transiently jammed clusters at the channel walls. 
In the early stage of the process, fast-growing hedgehog-like clusters are formed which are largely 
immobile. At later stages, most of these clusters dissolve and mobilize into nematized aggregates 
sliding past the walls. 

PACS numbers: 82.70.Dd, 05.40.-a, 61.20.Lc, 61.30.-v 



Swimming microorganisms, insects, birds and fish fre- 
quently move collectively in large groups with sponta- 
neous liquid crystalline order. Considerable recent re- 
search activity has been devoted to understand the origin 
of flocks and swarms in terms of simple models of self- 
propelled particles [j 0, 0, 0, 0, 0, 0, 0] . In these models, 
"active" rods are driven by their own motor along the 
rod orientation axis and dissipate energy in the suspend- 
ing medium. In fact, it has been shown that a local 
interaction between neighboring particles which favors 
mutual alignment induces a non-equilibrium phase tran- 
sition from a state of zero transport to another one in- 
volving a big swarm with cooperative transport [l] . This 
nematically ordered state has been analyzed with hydro- 
dynamic approaches Q and instabilities at large wave- 
lengths and giant density fluctuations were predicted 
[J, 0, HJ. In a model with explicit rod interactions, self- 
organized vortex clusters were found in which particles 
rotate around a common center 0, 0] . 

Recently self-propelled rods have been realized in flu- 
idized monolayers of macroscopic rods [9j and by means 
of vibrating asymmetric granular rods [1 01 ] . Giant den- 
sity fluctuation were confirmed @. In these realizations, 
the rod-rod interaction is the crucial input in contrast 
to collective motions of bacteria or low- Reynolds-number 
micro-swimmers [ll which are dominated by hydrody- 
namic interactions II, [Hj]. While most investigations 
have been performed in the bulk, mainly in two spatial 
dimensions (2D), there are few studies of self- propelled 
rods in confinement. Micro-swimmers moving in chan- 
nels have been investigated in Ref. [l4|. Self-propelled 
granular rods in circular-shaped cells were studied ex- 
perimentally in Ref. and clustering of particles near 
the cell boundaries was found [l(|. 

In this paper we study the effect of channel confine- 
ment on the collective behavior of self-propelled intcr- 



* Present address: Department of Chemical Engineering, Imperial 
College London, South Kensington Campus, London SW7 2AZ, 
United Kingdom 

^Electronic address: hlowen@thphy.uni-duesseldorf.de 



acting colloidal rods by means of Brownian dynamics 
computer simulation and dynamical density functional 
theory. Our motivation to do so is threefold: First, 
since every realization in nature is finite, there is a gen- 
eral need to clarify the collective behavior of active rods 
close to system boundaries, and, in particular, to fea- 
ture the transport properties of colloidal particles. Sec- 
ond, transport through narrow channels is omnipresent 
in many realizations, e.g. the propagation of microor- 
ganisms through veins and pores. Third, a predictive 
theory which starts on the level of the interactions and 
incorporates microscopic correlations is needed. 

A 2D slit will cause the rods to align along the chan- 
nel walls such that a nematized steady state with a local 
nematic director pointing along the wall direction is ex- 
pected for finite densities. However, for active particles, 
we observe a novel clustering phenomenon close to the 
walls occurring on top of a nematization. We find an 
aggregation of rods at the channel walls into transiently 
jammed clusters which exhibit a hedgehog-like structure. 
Most of these immobile "hedgehog" clusters show up in 
the early stages of the aggregation process. They can 
efficiently block particle transport through the channel. 
The clusters do not only form in 2D channels but also 
appear in thin capillary slits in 3D. For soft particle in- 
teractions considerable homeotropic order near the walls 
is found which is purely generated by non-equilibrium. 
The aggregation near walls and the stability of hedgehog 
clusters is also borne out by a dynamical density func- 
tional theory which we construct from the Smoluchowski 
equation HaHEE?!- 

The clustering process can be detected in experiments 
using catalytically driven nanorods [H, [l9[ and colloidal 
particles [20| in micro-channels which exhibit the Brow- 
nian dynamics of our model (2?| . A control over the 
hedgehog cluster formation has relevant implications for 
transport of active particles through confinements. 

In our model, we describe a 2D system of N charged 
Brownian rods of length Lq confined in a channel with 
a fixed width Ax = 5Lq. Rod- rod interactions are im- 
plemented via a segment model where each rod is parti- 
tioned into 13 equidistant segments [Hj]. Each two seg- 
ments from different rods interact with each other via a 
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FIG. 1: (Color online) (a-d): Simulation snapshots (cutout 
of the full simulation box) showing the aggregation process in 
a 2D system of self-driven rods with Fu = 5kBT/Lo, K-Lo 
If), (a) t = (isotropic), (b) t = O.lOrg (hedgehog), (c) 
t = 0.18i~b (nematized hedgehog) and (d) t = 0.73tb (wall 
aggregates). Clustered rods are shown in red, unclustered 
ones in black, typical hedgehog clusters are encircled. Arrows 
indicate the direction of Fu ■ (e-f): Wall clusters appearing in a 
3D slit capillary with width Ax = 5Lo and height Az — 0.3Lo 
at the same driving force. Shown are the projections onto the 
xy-plane. The parameters are kLq = 10, NLq/V = 5 (with V 



the system volume) , 
t = 0.15rs (f). 



t = 0.09r B (e) and kL = 5, NLt/V = 3, 



Yukawa potential [21] . The pair potential between two 
segmented rods at given center-of-mass positions {r,r'} 
and orientational unit vectors {u, u'} is then given by 

C/ rod (r-r ;u,u) = U Q }^}^ — , (1) 



with k the Debye screening constant and = |r — r' + 
(fjU — ^j-u')| the distance between segment i and j with 
—Lq/2 < £i/j < L /2. In all cases, U = 20 k B T with 
Boltzmann's constant and T temperature. The same 
Yukawa segment potential is used to model the rod-wall 
interaction: 



U waa {x,u) = U ^2 



exp[— KXi] 



(2) 



with Xi = \x + H.iU.\ and x the rod center-of-mass distance 
to the channel wall. This potential prevents any particle- 
wall overlap and enforces the rods to orient along the 
channel direction. 

To model an ensemble of self-propelled rods a constant 
force f || is applied along the main axis of each rod. The 
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FIG. 2: (Color online) Time-dependent profiles for the num- 
ber density p(x) (a) and nematic order parameter S(x) (b) 
for F|| = 10kBT/Lo, nLo = 10 obtained from simulation, 
(c) Adsorption Ap — f Q °[p(x,t) — p(x,0)]dx and (d) excess 

orientation AS = p{x)[S{x, t) - S{x,0)]dx/ dxp(x) 
showing the evolution of the rod structure with respect to the 
equilibrium initial state near the channel walls. 



system is kept at constant temperature T via the solvent 
and the dynamics is simulated by a finite-difference in- 
tegration of the overdamped Langevin equations using 
the short-time translational and rotational diffusion con- 
stants of a rod taking the limit of infinite aspect ratio 
[2l| . At any finite number density, the packing fraction 
of the system is then negligibly small and hydrodynamic 
interactions between the rods are of minor importance. 
The overall number density po — NLq/A (with A the sys- 
tem area) is fixed at po = 1 which corresponds to a 2D 
isotropic equilibrium bulk structure. Time t is measured 
in units of the Brownian relaxation time tb = Lq/Dj, of 
a single rod in terms of the translational diffusion coeffi- 
cient parallel to the rod d\,. All simulations were carried 
out with a time step St = 10 tb, with a periodically 
repeated finite system size of 125Lo along the channel 
comprising N = 625 rods. Rod clusters were identified 
using a distance criterion such that a rod pair is assigned 
to a cluster if their distance of closest approach is less 
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In all processes, the initial rod configuration is the 
isotropic equilibrium state for F\\ = 0. At t = 0, Fu 
is switched to a finite value with the forces pointing ran- 
domly along the main rod axis, as depicted in Fig.[TJi, and 
rods start accumulating at the walls forming tight clus- 
ters with a large fraction of rods pointing along the wall 
normal (Fig. [TJd). Owing to its symmetric, semicircular 
shape these wall clusters will be referred to as hedge- 
hogs. Since the net self-propulsion force is very small the 
hedgehogs are practically immobile which enables them 
to gather many other rods and grow rapidly. When the 
hedgehog has reached a certain critical size it becomes 
structurally unstable and slowly dissolves into a nema- 
tized cluster with nonzero cooperative motion (Fig. [it). 
Simultaneously, new hedgehogs may form at a vacant 
part of the wall. In the steady state (Fig. QJ1), the 
nematized hedgehogs have spread out along the walls 
and the rods move cooperatively in mono-layer aggre- 
gates displaying a characteristic tilt angle with respect 
to the wall and only few hedgehog clustering events are 
seen. Steady state density fluctuations are measured 
from a 2 = (A/ 2 ) - (M) 2 with M the number of rods 
in the left (L) or upper (U) half part of the simulation 
box. The ratios, <JL/o~L,eq = 11, o~T/o~u.eq — 20 (for 
Fu = Wh^T/La) confirm huge density fluctuations for 
confined driven systems. 

A simulation study in a 3D slit geometry (see Fig. QJ) 
shows that the formation of hedgehogs is robust against 
the rods moving out of the 2D plane. In Fig. [Tf we ob- 
serve that a reduction of the effective aspect ratio kLq 
leads to a transformation of the hedgehogs into linear 
homeotropic wall clusters consisting of arrays of parallel 
oriented rods. These type of clusters appear to be per- 
sistent in time, also in a 3D slit-geometry. Note that the 
homeotropic alignment is a purely non-equilibrium effect 
which is not imparted by the external wall potential. 

The behavior depicted in Fig. [T^-d also globally shows 
up in the profiles for the density p(x) and nematic order 
parameter S(x) across the slit, shown in Fig. [2] The 
latter is defined as 
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(3) 



with (-) x denoting a local canonical average at position 
x. S is close to unity if rods are preferentially parallel 
to the wall, whereas a negative value indicates a favor- 
able perpendicular alignment. The density profiles are 
obtained by a time-dependent canonical average of the 
microscopic density using about 500 independent pro- 
cesses. The sharp increase of the wall peak points to 
a strong aggregation of rods at the wall. The double 
peaked shape of p(x) at intermediate times can be at- 
tributed to the presence of both hedgehogs (associated 
with a sharp minimum in S(x)) and nematized wall clus- 
ters. Details of the aggregation process are depicted in 
Fig. [2t-d. There, the rod structure in the direct vicinity 
of the walls is shown in terms of the integrated differ- 
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FIG. 3: (Color online) (a) Average number fraction of clus- 
tered rods (N c /N) versus time. The solid curves correspond 
to an initial state of freely rotating rods, the points to an 
aligned state where Ui _L n for each particle i. (b) Time evo- 
lution of the hedgehog strength Sn(t) of a hedgehog nucleus 
(see Eq. ©) from dynamical density functional theory. 



ence of the profiles with respect to the equilibrium state 
at t = 0. The excess orientation (Fig. [2ji) varies non- 
monotonically with time and the minimum can roughly 
be identified with the point where the transient hedgehog 
clusters break up and transform into nematized ones (see 
Fig.[l|. From this we can then infer that: (1) the average 
lifetime of the hedgehogs becomes smaller upon increas- 
ing Fu and (2) the ones formed at smaller driving forces 
have an increased propensity to be perpendicular to the 
wall. A similar non-monotonic behavior is encountered 
in the average cluster size as a function of time shown in 
Fig. [3] For the largest forces, the maximum appearing 
in the cluster fraction roughly corresponds to the ones in 
Fig. [2}i which provides additional support for the tran- 
sient clustering scenario. 

We stress that at later times, say t > O.Itb, hedge- 
hogs do still appear (see Fig. [lb) albeit with a smaller 
probability. To ensure that the qualitative features of 
the clustering scenario do not dependent upon the initial 
structure, we also considered clustering starting from an 
equilibrated nematic system where all rods are enforced 
to be parallel to the wall. The result, included in Fig. [31 
shows that a similar clustering process, albeit less pro- 
nounced, takes place starting from a collection of nematic 
swimmers. 

Let us now turn to a microscopic theory for the aggre- 
gation process which is formulated in terms of the one- 
body density p(r, u, t) for a non-equilibrium ensemble of 
self-propelled rods. A general equation of motion for p 
is obtained by integrating the iV-particle Smoluchowski 
equation which leads to [13, H2, HH : 



dtp = V • D T • [Vp + pV (/W int + PU wm ) 

+ d r [d 2 ^ + d vP d v {pu int + 0U waU )] 



/3F|,pu] 



(4) 



where (3 = (k^T)^ 1 , V = {d x ,d y } and u = {sincp, costp} 
with ip the rod-wall angle. Furthermore, is the diag- 
onal translational diffusion tensor and Dr the rotational 
diffusion coefficient. The effective aspect ratio is fixed at 
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FIG. 4: (Color online) Same as Fig. [2] Results from dynami- 
cal density functional theory. 



kLq = 10. A connection with density functional theory 
is made by defining the non-equilibrium collective inter- 
action potential /3U- m t = Sf3F cxc [p)/Sp in terms of the 
equilibrium excess free energy functional F exc [23[. Here 
we use a simple Onsager functional 0, HU 



hit 



(r, u) = - J dr'du' (exp[-/3C/ rod ] - 1) p(r', u'), 



(5) 

which accounts for the rod interactions on the second- 
virial level. 

To describe the wall aggregation process we assume 
the density to be variable only along the wall normal n 
and the one-body densities for t > were numerically 
obtained from Eq. From Fig. 2] we see that the the- 
ory captures the initial stages of the aggregation process. 
The growth of the wall density peak signifies a strong 
adsorption at the wall while the concomitant decrease of 



S(x) points to the rods tilting away from the wall. 

To ascertain the stability of hedgehogs, we have solved 
Eq. ((3]) with full 2D spatial resolution focusing on the 
fate of a small semicircular hedgehog wall nucleus with 
radius R c at t — parametrized as follows: 



p(r,vj,0) 



' O eq (x) cx P[ a cos (f-7(r))] 



p eq (x,ip) 



r > Rc, 



(6) 



in terms of the local director angle 7(r) = arccos(j//r)+7r, 
the equilibrium density distribution p eq [with p eq (x) = 
J dipp eq (x 1 ip)] and local nematic order parameter a — 
100 {Io(ct) is a modified Bessel function). The hedgehog 
size can be measured in terms of the following weighted 
density 



K(t) 



drd^exp[-(r/i? c ) 2 ]p(r, ip, t), 



(7) 



with r the distance from the center of the hedgehog nu- 
cleus at the wall. Similarly, the hedgehog order follows 
from 



H(t) 



drdip cos(ip — 7(r))p(r, ip, i). 



(8) 



We may now define the hedgehog strength Su{t) = 
lZH(t) as a suitable parameter to assess the stability of 
the hedgehog nucleus over time. The result in Fig. [SJd 
shows that a small hedgehog of size R c = L is indeed 
stabilized provided that the self-propelling force is suffi- 
ciently large. Both the transient nature of the hedgehogs 
and the location of the cluster maximum as a function of 
the force are in agreement with the simulation results in 
Fig. 

In conclusion, we have shown that self-propelled rods 
confined in a channel exhibit an aggregation at the sys- 
tem boundaries. The aggregation process is dominated 
by a transient formation of virtually immobile clusters 
which possess a hedgehog structure. The results found 
in non-equilibrium Brownian dynamics computer simu- 
lations arc confirmed by a microscopic dynamical den- 
sity functional theory. Our predictions are in principle 
verifiable in experiments on colloidal suspensions [26j or 
catalytically driven nanorods [HI, [l!| . All these experi- 
mental systems are essentially realized in 2D. The aggre- 
gation phenomenon is important for blocking collective 
motion of active particles through narrow pores as en- 
countered in microfluidic devices. 
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